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ABSTRACT 


Experimental studies of isolated rotors in forward flight have indicated considerable 
potential for noise reduction through the application of higher harmonic pitch control. Tests 
to date also show that such pitch inputs can also generate substantial vibratory loads. This 
report summarizes the modification of the RotorCRAFT (Computation of Rotor 
Aerodynamics in Forward Flight) comprehensive analysis of isolated rotors to study the 
vibratory loading generated by representative high frequency pitch inputs. The RotorCRAFT 
code has been developed for use in the computation of such loading, and incorporates a 
highly refined rotor wake model to facilitate this task. This model discretizes the sheet of 
vorticity that trails from each blade by laying out vortex filaments along contours of constant 
sheet strength. In previous work this Constant Vorticity Contour (CVC) wake model has 
been found to be an important contributor to improved accuracy in the prediction of unsteady 
airloads. This report describes the extension of RotorCRAFT to include a variety of features 
needed for the study of vibratory response to higher harmonic pitch, blade/vortex interaction 
events, and the prediction of rotor noise. These include: application of arbitrary root pitch 
control; computation of internal blade stresses and hub loads; improved modeling of near 
wake unsteady effects; and preliminary implementation of a coupled prediction of rotor 
airloads and noise. Correlation studies are earned out with existing blade stress and 
vibratory hub load data. Sample calculations are executed to Illustrate the capabilities of the 
modified RotorCRAFT code and to point out possible directions for future work. 
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1. INTRODUCTION 


Efforts to predict and reduce the noise generated by helicopter rotors have been 
underway for several decades. Recent experiments (Ref. 1) have suggested that application 
of higher harmonic pitch control holds substantial promise as a technique for reducing the 
noise due to wake/blade interaction. However, results to date have also indicated that 
substantial vibratory loads can be generated as a consequence of such inputs, and it was 
judged important to understand the potential tradeoffs between possible noise reductions and 
increased vibration. 

To support the ongoing work in this area, the aerodynamic loading model in the 
RotorCRAFT (Computation of Rotor Aerodynamics in Forward flighT) free wake analysis 
of rotors in forward flight has been extensively modified and a wide range of demonstration 
calculations undertaken. Earlier work (Ref. 2) centered on the development of the Mod 0.0 
version of RotorCRAFT, a comprehensive analysis of rotor aerodynamic loading in forward 
flight using an advanced model of the vortex wake. The overall goal of the present effort 
was to enhance the capabilities of the baseline code to permit the computation of blade 
stresses and hub loads due to higher harmonic pitch inputs. As will be described below, this 
involved not only a substantial extension of the original structural model but also refinements 
in the modeling of unsteady loads generated by the near wake. Also, new output features 
have been added to permit the visualization of wake-induced velocity fields and the 
generation of unsteady airloads in a format suitable for input to NASA's WOPWOP code. 
The latter feature is present in a preliminary form only, however, and is suitable primarily for 
demonstration calculations of noise, not for data correlation efforts. A separate effort is 
presently underway to develop new, more efficient and accurate noise calculations. 

The RotorCRAFT code features unique capabilities that were judged to make a 
suitable platform for the research and correlation effort undertaken here. Among these is the 
Constant Vorticity Contour (CVC) wake model implemented in the code, a novel 
discretization of the wake into filamentary vortices representing constant-strength contours of 
vorticity in the vortex sheets that trail from each blade. The CVC model is noteworthy in that 
it eliminates the artificial distinction between vorticity generated by azimuthal and spanwise 
changes in the bound circulation on the blade (i.e., between what is conventionally termed 
"shed" and "trailed" vorticity) and in that it provides a visually meaningful representation of 
the wake. Previous work (Refs. 2-5) has indicated the success of this approach for obtaining 
improved predictions of wake-induced vibratory airloads relative to more simplified models 
using, for example, a single tip vortex. 

The computation of aerodynamic loads in the presence of the wake-induced flow field 
is carried out in the present analysis using a vortex lattice model of the blade. This model is 
appropriate for rotor calculations because of its potential for refined treatment of tip effects 
and wake/blade interaction. As will be described below, the development of the Mod 1.0 
variant has also involved studies of an improved model of unsteady near-wake aerodynamic 
loading, a feature that can be important for handling cases involving higher harmonic pitch 
control. Similarly, in cases involving high frequency pitch forcing, an appropriate treatment 
of the dynamics of the rotor blade is particularly important The present report documents the 
extension of the finite element model of structural deflection originally described in Reference 
2 to the problem of the prediction of internal blade stresses as well as forces and moments on 
the rotor hub. 

Finally, the results presented here include correlation studies of predicted rotor hub 
loads from the data available in Reference 1 as well as additional data from the literature on 
rotor blade stress computations. These correlations will be discussed after the presentation of 
information on the development of the present model in the next three sections. 



2. BACKGROUND ON THE ROTORCRAFT AERODYNAMIC MODEL 

A detailed discussion of the development of the baseline aerodynamic model is 
presented in Reference 2. This section briefly summarizes the major features of the technical 
background of the code. Before presenting this discussion, it is appropriate to briefly review 
past work on forward flight airloads in the literature. 

Scully (Ref. 6) developed a free wake analysis of the rotor using two free filaments (a 
rolled-up tip vortex and a diffuse inboard filament) to model the rotor wake. This distorted 
wake model succeeded in predicting certain features of the unsteady loading on rotors, 
though uncertainties about the proper modeling of close blade/wake interactions precluded 
accurate prediction of higher harmonic loads. Other broadly similar wake models nave been 
developed over the last fifteen years (Refs. 7-9). The seminal paper of Hooper (Ref 10) 
pointed out the inadequacy of these types of models in predicting the nearly impulsive 
airloading events that occur in a broad range of rotor designs, particularly m high speed 
flight. His results suggested that some aspects of the traditional approach to rotor wake 
aerodynamics must be inadequate, at least for the high speed case. Recent work by Miller 
and (Ref. 11) and Johnson (Ref. 12) has been directed toward repairing some of the 
shortcomings in previous single-tip-vortex models. Hooper’s results also motivated me 
development of the CVC model into a comprehensive free wake model suitable for 
application to rotors at both low and high forward speeds. The following sections describe 
the foundations of the CVC model and the methods used to incorporate it into the 
RotorCRAFT rotor analysis code. 

2.1 Constant Vorticity Contour (CVC) Rotor Wake Model 

The CVC method uses the curved vortex elements described in References 2, 13, and 
14 as the fundamental building blocks of a free wake model. The use of these elements 
contributes substantially to the accuracy and efficiency of the wake model and is particularly 
effective for representing the complex wake characteristics of forward flight. These 
complexities arise because the time-varying loading experienced by rotor blades generates 
both trailed and shed vorticity components in the wake. An idealized sequence of advancing 
side load distributions containing this behavior is sketched in Figure 2-1. A corresponding 
sequence of bound circulation distributions would exhibit similar behavior. Spanwise 
variations in the bound circulation are responsible for the "trailed" component of wake 
vorticity and temporal (i.e, azimuthal) variations are responsible for the "shed" component. 

Given a bound circulation distribution along the span, r(r,\|/), the wake that trails from the 
blade consists of a continuous sheet of vorticity whose strength as it leaves the blade is 

7(r,\|/). This can be represented by 

Y = Ysi + Ycj (2.1) 

where s and c are defined parallel to the span and to the chord of the blade, respectively. The 
intensity of each component of vorticity is related to the bound circulation by 

_ -anr,v) 

Ys " at (2.2) 


and 
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-dlXr.y) 
Yc= dr 


(2.3) 


while the magnitude of the vector is determined by 

Y = |y| = Vy s 2 + Yc 2 (2.4) 

Note that the strength of the wake sheet may be zero and that each component may reverse 
sign depending on the strength and rate of change (spatial or temporal) of the bound 
circulation. 

Figure 2-2 shows the idealized wake vorticity field corresponding to the load 
variations in Figure 2-1. The wake here is shown in terms of contours of constant sheet 

strength, i.e. lines of constant y are depicted (ignoring any vertical distortion). Clearly, only 
a finite number of contour lines are used in this figure, though a continuum of lines are in fact 
required to completely represent the wake sheet. Because increments in y are constant 
between each contour line, the amount of circulation contained between any two contour lmes 
on the sheet is constant. This increment in circulation is related to the bound circulation on 
the blade as follows: 


AT = I - y(r,y) dr = T(r+Ar,y) - T(r,\|0 

Jr (2-5) 


If the origination points of contour lines (denoted "release points") are spaced radially 
along the blade corresponding to fixed increments in bound circulation, then the radial 
distance between the contour lines will be a direct measure of the gradient in bound 

circulation along the blade. That is, if AT is held constant, the approximate relationship 


5f(r,\|/) g AT 

Ar(r,\jO (2.6) 

means that the contour spacing will be inversely proportional to the spanwise circulation 
gradient Thus, contour lines will be tightly spaced in the radial direction where gradients of 
bound circulation are high, as is typically true near the blade tip at most azimuth angles. As 
Figure 1 makes clear, however, realistic blade load variations do not always produce large 
bound circulation gradients near the blade tips. Gradients are often shallow at the advancing 
tip in high speed flight, and so the wake may take on a more sheet-like character, rather than 
rolling up into a concentrated tip vortex. Conversely, bound circulation gradients may not 
always be shallow in certain regions. Note in particular that the radial gradients near the 
blade root in forward flight may be relatively large compared to the weak gradients typically 
seen in the inboard loading distribution of hovering rotors. Thus, unlike the diffuse inboard 
wakes usually present on rotors in hover, concentrated vortical structures may be present 
trailing from the root region of rotors in high speed flight 
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Note that the tangent to the contour line at any point determines the direction of the 

local vorticity vector 7 , with the "trailed" and "shed" vorticity simply being components of 

the resultant vector 7 . While conceptually convenient to retain this terminology in some 
cases, it is desirable to employ a wake discretization that removes the artificial distinction 
between these two quantities and treats the vector vorticity field in a unified manner. The 
first step in constructing such a discretization is to recall that the contour lines can in fact be 
treated as vortex filaments. The natural step to take in representing the wake is thus to lay 
down curved vortex elements along the contours of constant wake strength, as shown in 
Figure 2-3. To conserve circulation in the wake, these filaments are assumed to contain all of 
the rolled-up vorticity trailed between adjacent release points. This means that all of the 

filaments will have equal strength AT, and the filaments have constant strength along their 
length. 


The shed and trailed vorticity is thus accounted for by the fact that the direction of the 
tangent vector to the resultant elements changes continuously even though the strength of 
each filament is constant. Using the CVC approach, on the average only half as many 
elements are required as in the straight-line lattice method, where the local vorticity vector is 
represented by two discrete straight elements. Thus, the use of curved elements means that 
fewer, relatively large elements can be used to accurately represent the vector vorticity field in 
the wake. These factors give the full-span CVC wake model a substantial advantage in 
efficiency relative to the lattice approach; since the model requires roughly half the number of 
elements, the number of operations required to compute the wake-on-wake velocity field will 
be reduced by approximately a factor of four. 

As noted above, by choosing each contour filament to be of equal strength the 
filaments are closely spaced in regions of high circulation, and are widely spaced in regions 
where the circulation is low. This relationship between filament spacing and sheet strength, 
and the fact that the filaments follow the actual vortex lines, means that the structure and 
relative strength of the wake can be seen visually. Since the computational filaments are 
aligned along the true vortex lines in the wake, this seems to be the most natural and accurate 
representation of the wake. Also, because the "shed" cross-filaments present in the lattice 
model are eliminated, an additional discretization along the wake length is not needed and the 
potential problems of handling these cross-filaments in a free wake representation are 
avoided. 

Once the wake is laid out, its spatial evolution proceeds in the same way as in 
traditional Lagrangian models. The key difference is that the wake now trails from the full 
span of the blade and is not a simple root/tip filament pair. Calculations discussed in 
Reference 3 have shown the important role played by the wake of the inboard portion of the 
blade. In addition, this full-span approach naturally allows for the appearance, evolution, 
and disappearance of the wake of opposite-sign circulation zones (e.g., near the tip) as they 
occur. Two new features of the wake representation that must be addressed are the formation 
of loops in the wake and the spanwise variation in filament release points as the blade rotates. 
Handling both of these features requires careful treatment of the numerical logic involved in 
filament generation and amalgamation. Current procedures in RotorCRAFT have proved 
robust in this respect. 

In its current form, the analysis starts with an estimate of the spanwise bound 
circulation distribution at a set of azimuthal locations, computed using a vortex lattice 
treatment of the blade loading in conjunction with a simple uniform inflow model. The 
spanwise release points for the individual discrete filaments are found by interpolation 
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between the fixed calculation points along the span. This must be done at each azimuthal 
station since the circulation distribution changes. Not only do the release points shift along 
the blade span but some of them also disappear and reappear as the maximum circulation 
decreases and increases. This disappearance and reappearance corresponds to the formation 
of loops in the wake. 

Although the use of the full-span CVC wake model is important to the computation of 
wake motion and wake-induced velocities in the region near the rotor, it is in general 
computationally inefficient to retain this model over the full length of a semi-infinite wake. 
Thus, provision has been made for using a set of successively simpler wake models that are 
invoked for regions far from the rotor disk. Reference 2 discusses the steps taken to allow 
efficient far wake modeling in detail. This reference and the RotorCRAFT Users Manual 
(Ref. 15) also describe the technical background of the present vortex core model as well as 
guidelines for its implementation. 

2.2 Flow Field Calculations 

As noted above, the Mod 0.0 version of RotorCRAFT was intended primarily as a 
tool to investigate unsteady airloading on rotors in forward flight. The role played by the 
wake-induced velocity field in such calculations is, of course, critical. RotorCRAFT (Mod 
0.0) allows for printing out the velocity induced at the control points of the vortex lattice, 
though it does not provide for evaluation of wake-induced flow fields at other points of 
evaluation. Such a capability is potentially very useful for analyzing important effects 
involving blade/vortex interactions, as well as for visualization of wake behavior. Thus, as 
described in Reference 15, users of RotorCRAFT (Mod 1.0) may select as many as 1200 
evaluation points off the rotor disk where the wake-induced velocity will be computed and 
printed out. These points are typically arranged in 'scan planes’, and the present 
implementation allows for as many as three separate scan planes containing 400 points. 

Figure 2-4 shows typical output produced using this capability. It illustrates two 
views of the instantaneous flow field in the vicinity of a UH-60 rotor in forward flight at 
advance ratio 0.3 and thrust coefficient 0.007. Figure 2-4a shows a top view of the rotor 
disk with a contour plot of the downwash contours in a plane normal to the rotor shaft and 
0.1R above the hub. (This particular rendering is in black and white, but much more 
instructive plots can be done in color with any number of current software packages). This 
plot shows many of the expected flow features associated with a typical rotor’s aerodynamic 
environment: a region of upwash near the front edge of the disk; a relatively large region of 
downwash in the first quadrant, downstream of the advancing side; steep velocity gradients 
near the edges of the disk, with the wake downstream of the advancing tip being more 
diffuse than that on the retreating side; and an area of relatively quiescent flow on the 
retreating side near the hub. Figure 2-4b shows a vector plot of the crossflow downstream 
of the same rotor, showing the expected swirling flow associated with the two tip vortices as 
well as a crossflow component on the inboard portion of the advancing side. These figures 
are merely illustrations of the capabilities present in RotorCRAFT (Mod 1.0); more detailed 
examinations of particular regions of the vortex wake can be found in References 16 and 17. 

2.3 Vortex Lattice Blade Model 

Many modem rotors feature complex planforms, and it is advantageous to use lifting 
surface analyses (a subset of aerodynamic panel methods) to analyze such designs. The 
particular model employed here focuses primarily on thin lifting surfaces with no side- or 
leading-edge separation (with the exception of tip effects in yawed flow discussed in Ref. 2). 
This is in the tradition of vortex lattice models first developed by Falkner (Ref. 18) and later 
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Figure 2-4. Evaluation of induced velocities in the vicinity of a UH-60 rotor at advance 
ratio 0.3 and thrust coefficient .007 . 




popularized by such researchers as Rubbert (Ref. 19) and Margason and Lamar (Ref. 20). 
The rotor blade is represented in the present work by a surface consisting of vortex 
quadrilaterals. Four constant strength straight-line vortices form the sides of each 
quadrilateral, except at the trailing edge of the blade where a modified vortex quadrilateral 
without the downstream line vortex is used. Also, as noted above, special modifications can 
be made to the side edge of the tip in cases where yawed flow causes part of the wake to be 
trailed from this edge. 


The lifting surface/vortex lattice analysis used here was adapted for use from the 
hover simulation described in Reference 21. This formulation allows substantial flexibility in 
the specification of the blade's planform so that complex designs may be accommodated. 
Currently, the lattice can be divided into as many as ten different regions, with separate linear 
distributions of twist, taper, and sweep within each. The spacing of the quadrilaterals in a 
vortex lattice analysis is an important consideration (Ref. 22). The judicious selection of the 
density, spacing, and orientation of the quadrilaterals can considerably enhance the efficiency 
and rate of convergence of the blade loading. 

The solution method used to find the bound circulation given this lattice is essentially 
a straightforward implementation of the classical approach described in the literature on lattice 
methods for fixed wing applications (e.g.. Reference 20). Given the location and orientation 
of each of the quadrilaterals on the blade, the velocity induced by the blade lattice on each of 
the control points is determined, assuming unit strength for each quadrilateral. Then the 
resulting velocity is resolved in the normal direction at each control point, yielding an array of 
influence coefficients relating the vector of bound circulations to the blade-induced 
downwash, w, at each control point: 


w 

where 


A 


Ag 



i , j = 1 Nq 


(2.7) 


( 2 . 8 ) 


Here, Nq is the number of vortex quadrilaterals on the blade. Note that the location of each 
of the quadrilaterals as well as the velocities above are defined in the blade reference system 
depicted in Figure 2-5. 


This array of coefficients is stored and the velocity induced at each control point by 
the free stream, the blade rotation and deflection (discussed below), and the rotor wake are 
summed and resolved into 'normal-wash' velocities that are then used to find the vector of 
bound circulation values on the blade that eliminates any flow normal to the blade surface at 
the specified control points: 

Ag+q.n = 0 ( 2.9) 

g = -A- 1 (q.n) ( 2 . 10 ) 

The bound circulation can then be used to solve for the forces and moments on each segment 
of the blade by applying the Joukowsky Law to each of the vortex quadrilaterals. This is 
summed to yield the total force on each blade in shaft coordinates: 
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Shaft Coordinates (X^Yj.Zs) 


Blade Coordinates (X,Y,Z) 



Figure 2-5. Definition of shaft coordinates and blade coordinates. 



Figure 2-6. Schematic of yawed trailers in the near wake of a rotor blade at 
\\f = 0° , advance ratio 0.4 . 
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F = (HI S + YJ S + TK S ) = SXF jk 

j k (2.11) 

Here, j is the summation index over the quadrilaterals on the blade, while k is the index of the 
blade edges. Moments exerted by the blade about the sectional quarter-chord reference axis 
can also be computed: 


M = (ll s + MJ S + NK S ) = £XFjkx?jk 

j k 


( 2 . 12 ) 


Moments about the flapping and lagging axes of the blade are taken to compute the gross 
blade motion, while pitching moments can be computed for each section to use for the 

calculation of torsional deformation. Here, rjk is the vector from the reference axis to the 

point of action of F. Tabulated information on pitching moment coefficients can also be 
used to determine the moments, as is described in Reference 2 . 


The forces predicted by the basic vortex lattice analysis must be modified to include 
predictions of viscous drag for performance calculations, as well as the effect of 
compressibility on the sectional lift and moment. Other modifications are necessary to 
account for stall, yawed flow, and unsteady effects in the near wake (Ref. 2). To introduce 
forces generated by profile (viscous and pressure) drag into the calculation, the only practical 
approach at present is the use of two-dimensional airfoil data. The form of data required for 
this coupling is tabulated values for sectional drag coefficients as a function of both Mach 
number and sectional angle of attack. 

Compressibility has an important effect on performance for Mach numbers typical of 
modem rotors. Some of this effect is captured by the inclusion of Mach number dependence 
in the 2D airfoil data used for profile drag coefficients. However, compressibility also has a 
significant impact on the lift generated by airfoils at specified angles of attack, and so its 
influence on thrust and induced power must be considered as well. For two-dimensional thin 
airfoils, treatments like the Prandtl-Glauert correction to lift curve slope are useful for 
including compressibility effects. In vortex lattice calculations, transformations that are 
similar in spirit but more elaborate in detail must be invoked. Currently, a transformation of 
the blade geometry is used that involves stretching the chord of the blade at a given radial 
station based on the local Mach number at that section (Ref. 21). The model also accounts 
approximately for tip relief effects due to airfoil thickness. 

2.4 High Frequency Unsteady Airloads 

As discussed in Reference 2, the original model of the near wake just downstream 
of each rotor blade consists of a set of straight vortex trailers that are extended back into the 
wake from the bound vortices on the blade. The orientation of these trailers may adjust 
with the local free stream direction at a given radial station (Fig. 2-6), thus providing a 
coarse approximation of the behavior of the near wake suitable for the resolution of low- 
frequency loads; however, any other structure of the shed vorticity within this overlap 
region is neglected. The effect of the shed vorticity in the wake far downstream of the rotor 
blade is of course captured with the full-span CVC free wake model. 

As part of the improvement of the aerodynamic model of the blade carried out for 
this effort, a more refined treatment of the unsteady loading due to shed vorticity in the near 
wake has been studied. The approach presently under investigation involves adjusting the 
downwash used to compute the loads on the vortex lattice with an approximate indicial 
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function that captures the near wake effects. In principle, the downwash could be adjusted 
by explicitly including vortex elements in the near wake, laying out components of shed 
vorticity parallel to the blade span to supplement the trailers shown in Figure 2-6. This 
approach involves a variety of computational difficulties, including the need to carefully 
position the discrete shed vorticity as a function of time step, and the necessity to incur a 
computational penalty associated with retaining high-resolution discrete vortex 
representations of the wake. 

Given this, it was judged desirable to investigate the use of an indicial function to 
compute the lift response of the rotor blade to the shed vorticity in the near wake. This 
approach is implemented in the vortex lattice method as part of the calculation of the 
induced velocity at a control point. Let w represent the downwash at a typical control point 
due to the influence of the wake and the blade motion. The modification required to take 
account of near wake unsteady effects can be formulated as a connection to the downwash 
w the resulting quantity w* will be used to compute the bound circulation distribution and 
can be stated as: 

w*(s) = w O rt (s,st,N) (2.13) 

where s is the nondimensional time and O rt is an indicial function that reduces to 1.0 in the 
case where near wake shed vorticity is omitted. The nondimensional time here is s=(Ut/c), 
where the reference velocity U is assumed to be 

U = fir + U„ siny (2.14) 

for a given radial position r and azimuth station V, and U«, is the free stream velocity. If 

\|/ ovl is the extent of the overlap near wake (Fig. 2-6), then sj is the nondimensional time 

required for the blade to traverse one azimuthal increment Ay. T is the time elapsed during 
one such blade azimuth increment and N is the number of azimuth increments contained 
within the overlap region. 

The indicial function assumes that the downwash at the control point ramps up 
linearly from zero over the time step T, reaching its full value at that point. The present 
vortex lattice scheme assumes essentially a step increase in downwash at each time step that 

does not generate any shed vorticity in the overlap near wake. The indicial function O rt is 
designed to add in the effect of this near wake while also including the effect of its 
truncation beyond y ovl so as not to double-count the wake that is captured by the CVC 
model. 


In the case of the present rotor wake analysis, the downwash is desired at discrete 
multiples of the time step T (which is assumed to be the ramp time St ) : 

w*(nst) = w O n (nst,ST,N) (2.15) 


The indicial function can be constructed from 
O rt (0,ST,N) = 0 
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and 


O rt (nsT,s T ,N) = <I>o(nsT,ST)+ X (ns T ,s T ,N) 

m=l 


(2.16) 


where 


O 0 (s,s T ) = 0 
for s < st and 

®o(s,ST)=l + £ln(&a±2) 

for s > St . The term d> m is determined by a recursion relation 


n-mN 

Orn (ns T ,ST,N) = X (^m-1 ([k+(m-l)N]s T ,ST,N) - d> m _i ([k-l+(m-l)N]s T ,ST,N)} * 
k=l 

3> w ([n-mN+l-k]sT,ST,N) 

(2.18) 


for n > raN and 

(ns T ,s T ,N) = 0 

for n < mN. The term <I> W in Equation 2. 1 8 is given by 


O w (s,st,N) = 0 
for s < st, while the term 


O w (s,St,N) = 


1 


(2N+1)s t L 


1 + 


-Lln(- S ' ST+2 ) 
st l s+2 Jj 


2 (2s + 2+ (2N-1 )s t ) (s - sr) 


(2N+1)s t (2s + 4+ (2N-1)s t ) (2 s + (2N-l)s T ) 
2 (In [(s - Sr+2) (2 s + 2+ (2N-l)sr)j -In [2(2N+l)s T ]) 

4(s+2) 2 + 4(s+2)(2N-1)s t + (2N-l) 2 s T 2 


+ 


(2.19) 


accounts for the wake truncation. 

The primary effect of this type of indicial function is to introduce a phase delay into 
the response of the aerodynamic loading while also decreasing the magnitude of the response 
to a given upwash or blade motion input. This method has undergone preliminary testing 
and has so far produced qualitatively reasonable results. However, it is has not yet been 
thoroughly validated and so has not been used in the calculations described below. The 
present version of RotorCRAFT (Mod 1.0) contains provisions for activating the near wake 
model when the final version becomes available. 
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3. BASELINE MODEL OF BLADE DEFORMATION AND MODAL PROPERTIES 

The structural properties of the helicopter blade are clearly critical to the evaluation 
of the response and performance of a helicopter in forward flight. This is particularly true 
of the study of rotor vibratory hub loads undertaken here. Thus, one of the major efforts 
within the development of the RotorCRAFT code has been to incorporate a realistic finite 
element (F.E.) representation of the blade. The primary emphasis of the Mod 0.0 version 
of RotorCRAFT was to predict unsteady aerodynamic loading that leads to rotor vibration 
The pertinent final report. Reference 2, contains a full discussion of the major features of 
the finite element model and details of the formulation of the dynamic modeling and its 
capabilities and limitations. As noted above, this present effort involved the 
implementation of computations of hub loads and blade stresses. Thus, m the discussion 
that follows the major features of the finite element model shall be concisely reviewed while 
the enhancements implemented in RotorCRAFT (Mod 1.0) shall be discussed in greater 
detail The description of its coupling with the aerodynamic model discussed above is 
reserved for Section 4. Further description of the inputs required for this portion of the 
analysis are given in Reference 15. 

3. 1 Summary of the Finite Element Structural Model of the Helicopter Blade 

A fourteen degree-of-ffeedom finite element (F.E.) model is used to represent the 
helicopter blade extension, twist and transverse bending displacements. Stiffness and mass 
properties for each element are computed from the cross-section geometry and material 
properties supplied by the user. Blade rotation forces and geometric stiffening effects are 
accounted for in the construction of the elements. The resulting element matrices are then 
assembled and any constrained d.o.f. eliminated to finally yield global mass and stiffness 
matrices for the complete blade structure. The approach taken is similar to previous 
implementations of F.E. methods for rotorcraft applications, such as Reference 23. The 
resulting mass and stiffness matrices serve as inputs to a generalized eigenvalue problem 
which is solved by a Jacobi iterative technique to obtain the modal frequencies and shapes 
for the F.E. model. The modal properties are dependent upon the frequency of blade 
rotation since the geometric stiffening is proportional to the square of this frequency. 

The transfer of information between the structural and the aerodynamic models in 
RotorCRAFT occurs chiefly via the modal properties of the blade. The mode shapes are 
used to compute generalized modal forces from the distributed aerodynamic forces, and 
these modal forces drive the corresponding modal responses. These responses in turn are 
used in conjunction with the mode shapes to determine the instantaneous displacements and 
velocities at any point along the blade. Finally, the blade deformations and deformation 
rates provide the necessary reformation to update the flow field and aerodynamic forces, as 
described in the previous section. Hence, one iterates toward a steady state solution where 
the blade displacements and their rates are constant at each azimuthal location. The 
remainder of this section contains a brief review of the main features of the F.E. model 
with greater emphasis placed upon the computation of hub loads and blade stresses. 

3.1.1 Blade Geometry 

Each of the blade segments defined in the RotorCRAFT blade geometry input is 
subdivided into structural finite elements in the manner specified by the user. The global 
axes for the assembled blade are denoted by XYZ in Figure 3-1, which shows the 
relationship between shaft, wind, and blade axes. Local axes, xyz , are defined for each 
element such that the x-axis coincides with the elastic axis, or the line of shear centers, of 
the element Axes y and z lie in the cross-sectional plane of the element and are related to 
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a) Rotation of wind axis system about Y w to reach shaft axis system 
(Z $ parallel to rotor shaft) 



Z s > Z 

b) Rotation of shaft axes about Z s to reach blade axis system (X,Y,Z) 
Figure 3-1. Axis system definitions. 



the global Y and Z axes via a local to global transformation matrix, [T rot ]. This matrix is 
constructed from blade segment anhedral, sweep and collective specifications supplied by 
the user in the blade geometry input file as is -explained in Reference 2. The glob 
coordinates, XYZ, of any point on a given F.E. with local coordinates xyz are then given 

by: 




[Trot] 


m 


(3.1) 


where the coordinates, XoYoZq , are the global coordinates of the origin of the local axes 
here taken to lie on the elastic axis at the end of the element nearest the rotor hub. The same 
transformation matrix is used below in referring local deformations to the global axes 
Note that the rotations due to deformation can also be referred to the global axes using this 
transformation since the deformations are assumed small and the associated rotations thus 
commute. 

3.1.2 Element Degrees of Freedom 

The specification of the shape functions and the fourteen degrees of freedom of 
each element is summarized here. Each element has two end nodes and one node at its 
midpoint, as shown in Figure 3-2. The degrees of freedom correspond to translational and 
rotational deformations at these nodes. The deformation of the element at any point is 
estimated by interpolation of the nodal displacements using the shape functions. Let u,v 

and w denote the displacements along the local x,y, and z axes respectively and 0 the 
twist deformation about the x axis. Then the generalized displacement vector is defined as: 


Lag: 


Flap: 


Twist: 


Axial: 



(3.2) 


where the subscript (»), x denotes the derivative with respect to the local x axis coordinate 

and A is the element length. Thus qi and q2 refer to the displacement and 
corresponding slope due to bending in the y-direction at the left hand node. The 
corresponding right hand node deformations are q 3 and q 4 , and so forth for the other 
displacements. Note that the slopes, v, x and w, x , can be regarded as a small positive 
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rotation about the local z-axis and a small negative rotation about the local y-axis, 
respectively. 

The transverse displacements, v and w are interpolated using cubic Heimitian 

polynomials and quadratic polynomials are used to interpolate the torsional and axial 
deformations. This is the simplest element interpolation scheme yielding a consistent 
formulation for coupled torsion-bending (Ref. 23). Specifically. 


v(x) = {O 3 } 


t 


qi 1 
q2 , 

Q3 

q4 J 


w(x) = { 03 } T< q7 

. qs - 


0(X) = {O 2 } T { qio 1 

l qn J 
rqi 2 

u ( x )={< D 2 } T <113 

l qi4 


(3.3) 


where, 


{<h 3 } = 


< 




1 - 3^ 2 + 2£ 3 
($-2$ 2 + ^ 3 )A 
3^ 2 - 2^ 3 
(- ^ 2 + ) A 


"V 

' 1 - 35 + 25 2 ' 

> ; {02}=^ 

45 - 45 2 • 


-5 + 25 2 J 


(3.4) 


and £=x/A . The preceding relations may be expressed in compact form as: 


v(x) ' 
w(x) I 

0(x) 

u(x) j 


[<&] (q) 


(3.5) 


where {q> is the vector of generalized displacements and [ ] is a (4 x 14) matrix 

appropriately constructed from Equations 3.3 and 3.4. 

The global degrees of freedom are obtained by resolving the deformations along the 
global axes using the transformation matrix derived previously. At an end node, all of the 
three translational and three rotational d.o.f. are available (since the slopes of the transverse 
bending displacements correspond to rotations). Thus, the translation between the local 
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element d.o.f.'s, {qj^, . and die global ones, {q} G , is achieved using a straightforward 
transformation involving the matrix, Trot • 


IS) =[Troti { l } ; { I } =tTrot, { } &6) 

where Rx » R Y > and , are rotations due to deformation about the global XYZ 

axes respectively. At the mid-node the global deformations are defined to coincide with the 
respective local ones, i.e., 


fqio 1 = fqi o 1 

lqi3i G l<ii3j L 


(3.7) 


This both simplifies the transformation and results in an orthogonal element transformation 
matrix, i.e., if the elemental transformation which will be composed of elements of [T rot ] 

is denoted by [T GL ] so that {q} G = [T GL ] {q) L then [Tol]' 1 = [T GL ] T . 

3.1.3 Derivation of the Element Strains and Stresses 


In order to compute the elemental stiffness matrix the strains arising from the 
preceding displacements must be evaluated. The nonlinear expressions for the strains are 
stated: 


e xx = u, x + ( ^Q, x ), x + \ ( ti 2 + c 2 ) e , x 2 

+ \ v *x 2 - v, xx { T] cos(P+0) - C sin(P+0) } 

1 2 

+ 2 w ’ x " w, xx ( C cos(p+0) + T] sin(P+0) } (3.8a) 


Txtj - 2 e xr| = ( 'i'.ri - C ) ( 6.x + ®nl ) (3.8b) 

YxC = 2e x £ = (^ + Ti)(0, x + 0 nl ) (3. 8c) 

and all other strains are assumed zero. Here, T) and ^ , are the principal axes of the cross- 
section which are oriented such that: 

f"nC dA = 0 (3.9) 

In RotorCRAFT, the angle p between the local y and T) axes is specified by the user in 
the cross-section input file. In Equations 3.8, 0 n j is a nonlinear 2 nd order torsion term, 
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and ¥(x,ti,C) is the Saint Venant warping function expressing the out-of-plane 
displacement, u wai p , due to torsion: 

u waip (x,ri,0 = 'F(x,n,Q e, x (3 - 10) 

The linear expressions are easily obtained from above. 


fc xx 


Txri 


= u, x + ( 'PO.x )»x 

- v, xx { T) cosP - C sinP } - W, xx { C cosp + i) sinP } 

= (*V; )e, x ; Yxc = (^<+n)0,x 


(3.11a) 

(3.11b,c) 


These strains are expressed in terms of the vector of generalized d.o.f., {q} , apd sh ^P e 
functions and their derivatives w.r.t. x , by substituting for the occurrences of u , v , w , 

and 0 using Equations 3.3-3.5. This results in: 


b XX 


f - [ TjcosP - CsinP ] {$3 } ^ 
- [ Ccosp + Tjsinp ] {O3 } 

¥, x {02> +^{^ 2 "} 
{® 2 ) 


(q) 


(3.12a) 


Yxti 
YxC 






(3.12b) 


The corresponding stresses are derived from the Hooke s Law. 
“ = E £ xx 


} xx 


o xT1 = G Yxri ; °xC = G YxC 


(3.13a) 

(3.13b,c) 


These expressions are used both in the construction of the stiffness matrices and in the 
evaluation of the stresses at specified locations on the blade. The former involves integrated 
quantities of the shape functions and cross-section and material properties, wh ^ s , t ^ e ^ t 
requires knowledge of the local coordinates where the stress is to be evaluated and warping 
function related quantities evaluated at the required x-locauon. 

3.2 Evaluation of the Stiffness and Mass Matrices 

The discussion above defines the relationships for stresses and strains for the current 
F E formulation. Reference 2 describes the procedure for obtaining die element stiffness 
and mass matrices via application of the Extended Hamilton s Principle, and the assembly 
of the corresponding global matrices. Briefly, the process begins by deriving expressions 
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for the kinetic and strain energies of each element Any remaining forces not accounted for 
in the kinetic and strain energy expression, in particular the forces and geometric stiffening 
due to blade rotation, are represented by an external virtual work term. Variational calculus 
is then applied to the Lagrangian, £ = (Kinetic Energy) -(Strain Energy), which yields the 
equations of motion that govern the blade dynamics. Since in this case die variation of the 
strain energy is equivalent to the internal virtual work which can be written down directly, 

W 1 = J <*xx * 8e xx + o xT1 • 5y XT j + a x £ • 5y x ^ dV (3.14) 

the preliminary step of obtaining an expression for the strain energy is circumvented. The 
mass matrix is derived from those components of the kinetic energy expression containing 
products of the deformation velocities. The stiffness matrix is more complicated containing 
contributions from the strain energy, kinetic energy and from the external virtual work term 
representing the geometric stiffening. The volume integral is conveniently represented by 
an integral over the blade cross-section area, followed by an integral along the element 
length. Performing the area integration produces cross-section parameters such as area, 
first and second moments of area, torsion constant, etc., which in RotorCRAFT are 
supplied by the user in the blade cross-section input file. The integration along the element 
length is accomplished by 4 point Gaussian quadrature. 

3.2.1 Assembly of the Global Mass and Stiffness Matrices 

The element matrices obtained above are assembled into corresponding global 
matrices for the complete helicopter blade. The assembly process involves two sub- 
procedures: the first involves referring the elemental matrices and nodal forces to the global 
axes, and the second defines the array indexing that relates the local degrees of freedom for 
each element to the global ones. 

Rotation of the element matrices and nodal forces into a global coordinate frame is 
accomplished in the standard manner: 


[ K global ] = ITgl] [ K local ] [TglJ T 
[ Mgiobal ] = [TGLi [ M local ] [T G l] T 

{^global = [TGLi (FJlocal (3.15) 


as may be easily verified by noting that the kinetic and potential energies and the virtual 
work are independent of the choice of reference frame. Here, [TqiJ is the transformation 
matrix described in Section 3.1 relating local and global generalized d.o.f.: 
{ ^global) = [ T GLHqiocal) • 

The blade elements are then laid end to end in sequence from blade root to blade tip. 
Global deformations are defined as outlined in Equation 3.6. However, the ordering of the 
global degrees of freedom is different from the local ones. Each element degree of freedom 
is associated with a global one via an indexing array or splay matrix, C(k,ie) , where k is 
the local degree of freedom (k=l,2,...14) , ie is the element number, and C(k,ie) is the 
global degree of freedom. Having specified a C(k,ie) for each element then the global 
matrices may be constructed by ‘splaying’ components of the elemental mass and stiffness 
matrices into their corresponding positions in the global matrices. For example, the [i,j] 
entry of the elemental stiffness matrix for finite element, ie, is added to the 
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r ca ie) C(i ie) 1 entry of the global stiffness matrix. In like manner, the global nodal 
force vector is built upffom nodi forces for each element Section 5 of Reference 2 details 
the definitions of the local and global degrees of freedom. 

To this point in the development, the boundary conditions at the root have not been 
imposed, and some of the root degrees of freedom are eliminated when these are applied. 
For articulated blades, it is implicitly assumed in RotorCRAFT that die blade is freely 
hinged in both flap and lag directions, but that the remaining degrees of freedom at the root 
- the three translational displacements and the twist about the X axis - are constrained. 
Thus, the rows and columns of the global mass and stiffness matrices corresponding to 
these four degrees of freedom are removed prior to conducting the modal analysis. For 
cantilevered blades, by definition all root deformations are zero and thus all six degrees of 
freedom at the root are similarly eliminated. However, although the rows and columns 
corresponding to the constrained degrees of freedom are removed when computing the 
natural blade modes, they are employed when determining the inertia forces due to cyclic 
pitching, and also when calculating hub loads. 

Let the global mass and stiffness matrices for the blade be given by: 



Moc • Mcr 


Kcc ' Kcr 

M = 



; K = 

... . . ♦ « • • 


i 

* 

i= 

i 


r 

■r 


(3.16) 


where the matrices have been partitioned into submatrices associated with constrained and 

unconstrained nodal degrees of freedom respectively. Submatnces, [Mj-f] and |RrrJ are 

the mass and stiffness matrices for unconstrained deformations and are used in computing 
modal properties. The constrained degrees of freedom are associated with the translational 
deformations, U, V , and W, and also the twist rotation, R X , at the blade root. For 
cantilevered beams, the constrained degrees of freedom also include the rotations, Ry and 
Rz . Of particular interest are the matrix partitions associated with the root twist 

displacement, 0 C , which shall be required for computing the effects of cyclic pitch. 
Although the twist deformation is zero at the root, it is clear that there is a finite rigid Ixxly 
rotation occurring at the root due to cyclic pitching and this will introduce an associated 
inertia tom into the equations of motion governing the remaining degrees of freedom, i ne 
equations of motion obtained from the Extended Hamilton’s Principle are summarized here 
(see Ref. 2 for details) as: 


[M*.: 


r ^ 

• • 

e c 

0 c 

... 

+ ... 

•• 

{ qr J 

l q r J 



1 = {F aeroj + {protj 


(3.17) 


where [M^ c ] and [K^J are the column of the mass and stiffness matrices associated with 

the root twist {F 3 " 0 } are the nodal forces due to aerodynamic loads and {F rot } the nodal 
forces associated with the blade rotation. This last term is time invariant and therefore need 
only be computed once for a given blade rotation speed. Note that none of the other 


22 


constrained degrees of freedom, U, V, W, (and Ry , R z , for cantilevered blades) are zero 
and so do not appear in Equation 3. 17. The only non-zero root displacement is 0 C . 

The generalized displacement vector, {q,.}, includes both deformation effects and 
the rigid body rotation due to cyclic pitch. It is convenient to define deformation variables 
{q} by subtracting the rigid body motion arising from cyclic pitch. These represent 
deformations relative to a set of axes rotating with the blade root. 

{q} = (qr) - { r 9c> e c ; (q) = (q r > - {re c >e c (3.18) 

where {rQ c } is the vector of generalized displacements arising from unit root pitch. The 
components of {tq c } may be computed explicitly from the geometry as the global 
displacements and rotations due to unit cyclic pitch: 

U = 0 ; V = -r z ; W = r Y 

e x = l; R y =0 ; R z = 0 (3.19) 

Here r Y and r z are the global Y and Z coordinates of the elastic axis at a given location 
along the blade. One can also compute {rQ c ) from 

( r © c ) = - IK,,]- 1 [K^J (3.20) 

which is obtained simply by noting that rigid body rotations do not contribute to material 
strains. Strictly speaking, the latter expression is an approximation valid for small pitch 
angles. The effects of geometric stiffening are in fact dependent upon the orientation of the 

structure. Thus, changing 0 C alters the geometric stiffening entries of [K^]. Large cyclic 
inputs imply that the rotating frame acquires significant angular velocity and acceleration 
components along the global X-axis in addition to the steady angular velocity arising from 
blade rotation. A full analysis accounting for the various cross-products of translational 
and arigular velocities would be necessary, as is done say, in slewing manoeuvres of 
spinning satellites. In the present case, such an analysis is unwarranted in light of the 
linearization assumptions already made in the F.E. analysis and the small cyclic pitch 
amplitudes under consideration here. 

Rearranging Equation 3.17: 

[M*](q'r)+[KJ<qr)= (F 30 ™) + {F 101 } - ( [M rtc ] e c + [K^.] 0 C ) (3.21) 

and distinguishing between the blade deformations and rigid body displacements via 
Equation 3.18 

[M rr ](q)+[K 1T ](q) = (F 26 ™) + (F 101 ) - ( [M^+CMJIr*) )S C 

- ([K rtc ]+[K„](re c ))0 c 

= (F"°) + IF"*) - ( [Mrffc] + [Mn-Kr^J )S C (3.22) 
where use is made of Equation 3.20 in the final equality. 
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3.2.2 Computation of Modal Properties 

The dynamic analysis in the RotorCRAFT code is formulated in terms of the 
response modes of the blade. Thus, the information contained in the blade mass and 
stiffness matri ces is used to obtain a set of mode shapes and their corresponding modal 
frequencies and masses. The modes are obtained by solving the generalized eigenvalue 

problem: 

[Kn-] { xi } = ©i 2 ^] { xj } (3,23) 

where C0j is the natural frequency of the i-th mode, and {xj} is the corresponding 
eigenvector, or, mode shape. 

Equation 3.23 is solved using the an efficient and accurate generalized Jacobi 
iteration method described in Reference 24. The modes are then ordered in ascending 
frequency, and seated into the dominant type of deformation present in each mode: bending 
(flap and lag), torsion or extension. The number of modes of each type retained m the 
dynamic analysis in RotorCRAFT is specified by the user. 

Having obtained a set of modal frequencies and shapes, the generalized 
deformation, {q}, may be expressed in terms of modal amplitudes, {t)}, as: 

{q} = [X](n) (3 - 24) 

where [X] = [ (x;} ] is the set of column eigenvectors of the retained modes. Equation 
3.22 may then be transformed into the modal variables by pre-multiplying by [X] T and 
substituting {q}=[X]{r|}: 

[X] T [ M n ][X]{n } + [X] T [ K„ ][X]{ti} 

[X] T [ {F 3610 } + {F"*} - ( [MrfJ +[M IT ]{r 0c } ) 0 C ] 


or, GM(i) [ rii + £»i 2 Tli ] = F AER0 (i) + F R0T (i) - MRC(i) 0 C 

(3.25) 

where, 


GM(i) = {x i } T [M rT ]{x i ) 

(3.26a) 

°H - GM(i) 

(3.26b) 

F AER0 (i)={x i } T {F 3 * 0 } 

(3.26c) 

F R0T (i) ={x i } T {F rot ) 

(3.26d) 

and, MRC(i) = {^([M^] + [M^Hr^]) 

(3.26e) 
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The generalized mass, GM(i), is determined by the normalization method of {xj}. In 
RotorCRAFT, the modeshapes are normalized so as to have tip amplitude equal to the 
blade radius. Thus, for the first flap and lag modes, GM(i) is approximately the second 

moment of inertia of the blade about the hinge and the modal amplitudes, {t| } , correspond 
to flap and lag angles. 

3.3 Computation of Loads and Blade Stresses 
3.3.1 Hub Loads 


Hub loads are computed 
constrained nodes: 


from the equations of motion corresponding to the 


{F c } 


[Mce^Mcrl 


kj 


0 


r + 


q r 


c 

Or 


(3.27) 


where [M C0C ] and [K c q c ] are the column vectors of [M cc ] and [K^] respectively 
corresponding to the root twist deformation. The terms: [K CT ], [M cr ], [M^J, and [K c0c ] 

are directly extracted from the assembled finite element matrices prior to the constraint 
deletion procedure. The constraint forces, (F c ), are interpreted as the forces that must be 
imposed upon the blade in order to effect the constraint conditions. Therefore, the equal 
and opposite hub forces (experienced by the hub) are {F hub } = - (F c ). If qj-(t) is known, 
then one can in principle compute the constraint forces at the blade root. From the same 
arguments following Equation 3.20, one may rewrite the last term involving the stiffness 
matrix to obtain: 


{F c } =[Mc0 c :Mc r r 
with {q} defined in Equation 3.18 


6 


q r 


r + [KcrHq) 


(3.28) 


Potential difficulties in the computation of hub loads may arise when higher modes 
modes are neglected and q^t) is constructed from a reduced set of modal responses. 
Neglected modes which are unimportant when predicting aerodynamic loads and blade 
deflections may be quite significant in hub force computation. Consider for example the 
extension response of the blade. In most rotor blades, the frequency of extensional 
vibration is an order of magnitude greater than the torsional, flap and lag frequencies. 
Furthermore, the aerodynamic force component along the blade is usually negligible and 
the contribution of the blade stretching dynamics to the overall motion can thus be ignored. 
However, the forces due to blade rotation act primarily in the direction along the blade and 
thus the chief hub load arising from blade rotation acts in the X-direction which is precisely 
the axis along which the extensional motion occurs. So if extension modes are neglected in 

the analysis and hub forces are computed from {q}=[X] {r| } and Equation 3.28 then the 
dominant hub force in the X-direction due to blade rotation will be completely neglected. 
The fundamental issue is that the contributing forces may have significant components that 
are orthogonal to the space spanned by the retained modal set 
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To remedy this situation, it has been found useful to compute certain of the hub 
force contributions directly, and to then superimpose the effects due to the modal 
accelerations. The simplest derivation is to suppose one has retained all of the modes in the 
analysis so that the original FJE. model and the modal analysis yield identical results, to 
then evaluate all those forcing terms that can be computed directly, and to finally 
approximate the remaining terms with a restricted mode set. 

Rewriting Equation 3.28 

{F c } = [M C0C ]0 C + [Merit q r } + [^Hq) (3.29) 

and substituting for {q} and { q r } using Equations 3.18 and 3.24, 

(F c ) = [Mer][X] {fj} + ( [M C0C ] + [M^] {r 0( .} )e c + [K^HX] {q } 

(3.30) 

where {r 0( .} is computed from Equation 3.20 and [X] follows from the eigenanalysis, 
Equation 3.23. Using the modal equations, Equation 3.25, to solve for {q}: 


{q} = [Dl'^tF^ 0 } + {F rot } - {MRC} 6 C ] - [diag{ 1/C0j 2 }] {q} 

(3.31) 

where [D] = [ diag{ CDj 2 GM(i) ) ]. Substituting for {q} in Equation 3.30 using Equation 
3.31 and collecting terms: 

{F c } = ([Mcr][X] - [Kcr][X][ diag{ 1/(0? }]){q} 

+ ( [M^] + [McrHr^} - pyiXHPr 1 {MRC})e c 

+ [Ker] [X] [D] * 1 [ { F^n 0 } + {F R0T } ] (3.32) 

The last two terms: 

[K ct ][X][D]' 1 [F AERO (i) + F R0T (i)] = (K cr ][X][D]- I (X] T [(F a "») + (F rol ) ] 
but since 

[XltDl'^X] 7 = [X]([X] T [K rr ][X])’ 1 [X] T = [Krr]- 1 (3.33) 

Equation 3.32 becomes, 

(F c ) = ([Mc r KX] - [KcKXlt diag( l/ 0 ij 2 ]])(n ) 

+ ( [Med,] + [Mcrl(ree) - [KeJPqtD]- 1 [MRC] )0 C 

+ [Ke r ][K rr r 1 [(F“ ro ) + (F™) ] 
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where 

[MFC]{r[ } + {MCYC}9 c + {F^ 00 } + {F ROTC } 

(3.34) 

[MFC] 

= [MJfX] - [KJtXlfdiagf 1/cDi 2 }] 

(3.35a) 

{MCYC} 

= [McOcl + - [KjrapD]- 1 {MRC} 

(3.35b) 

jpAEROC 

} = [KJtKJ- 1 ^ 0 } 

(3.35c) 

{F rotc } = [K cr ][K rr ]- 1 {F rot } 

and {MRC} is the vector with components defined in Equation 3.26e. 

(3.35d) 


Fortunately, the forces due to steady blade rotation are constant, and thus {F R ^ T ^} 
is constant and need only be computed once. Furthermore, it does not contribute to the 

modal accelerations but in effect induces a steady state bias in the {t|}. This implies that no 
error in the rotational force contribution to the hub forces is incurred when using a reduced 

set of modes since the remaining quantities in Equation 3.34, {if }, 0 C and {F aero } are 
independent of the rotational force. In other words, all of the rotation force contribution is 
contained in {F R0TC }. 

The term. 


nc c j[K„r 1 =([K, r r , (K rc j) T 


(3.36) 


can be shown to be the transpose of the array whose k-th row corresponds to a unit 
displacement in the negative k-th Cartesian direction. For example, if the beam were given 
a unit displacement in the negative Y direction (corresponding to k=2) then the resulting 

vector (q}= ( [KJ^IX^.] ) k=2 . This is completely analogous to the definition of {tq } in 
Equation 3.20. Here, Cartesian directions, k= 1,2,3 correspond to displacements along the 
rotating blade X,Y,Z axes respectively, and k=4,5,6 correspond to rotations about those 

axes. Therefore, {pAEROCj - s g[ m pjy summation of the aerodynamic nodal forces in 
the k-th Cartesian direction which is obtained in RotorCRAFT by summing the 
aerodynamic forces and moments about the hinge directly. So, 


tip 

{ pAEROC j = J (pjaero ds 
cutout 



(3.37) 
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which is equivalent to the aerodynamic force and moment vectors defined in Equations 
2.1 1 and 2?12 when referred to the rotating Made. axes. Equation 3.34 « inAe same 
form when a reduced set of modes is employed. {MCYC} is evaluated in teims of the 
reduced mode set rather than the complete set since the cyclic pitch acceleration influences 

the modal accelerations so that the combination of the terms involving {rf } and 0 c m 
Equation 3.34 is consistent. Computing the aerodynamic loading via Equation 3.3 >7 'rather 
than in terms of the reduced modes is consistent smce the error incurred is of the same 
order as that due to neglecting higher modes in the first place. Furthermore, the hub loads 
averaged over one blade rotation equal the correspondingly averaged aerodynamic forces, 
as is expected theoretically. This is not generally true if the hub loads are computed based 

on the reduced set of aerodynamic modal forces, {F^ 0 }. The computed hub loads are 
seen to correspond to a blade subject to the true aerodynamic loads but with the 
accelerations determined from the restricted set of modes. 

3.3.2 Blade Stress Calculations 

The linearized longitudinal and cross-section shear stresses follow directly from 
Equations 3.12 and 3.13: 


J XX 


= E 



= G 


- [ T\cosP - CsinP ] {O 3 } 

- [ £cosp + T|sinP KO 3 } 

y.x {^ 2 } +^{^2"} 

{<*> 2 ) 

T 

0 
0 

{<*2} 

0 


{q>L 


(3.38) 



(q)L 


(3.39) 


where the subscript, L, denotes local generalized coordinates. Both the longitudinal and 
shear stresses are functions of the area coordinates, t\ and C, which must be specified. 
Furthermore, the stresses are most conveniently referred to the local reference which 
is desirable since the remaining components of the sttess tensor arc zero by as • 
Referring stresses to the global blade fixed axes, though straightforward in pnnciple, 
results in an awkward description of the stress distribution where, in general, all six of the 
stress components are non-zero. 

Note that while this treatment is adequate for straight, uniform blades, 
inconsistencies in the stresses referred to the local element axes arise whe u there is a 
discontinuity in sweep, anhedral or other structural or geometrical property between 

elements, i.e., a plot of, say, Oxx as function of x for a general blade topology ^1 exhibit 
discontinuities. In theory, near such points, the stress tensor is full and the simplified beam 
theory is not valid. However, by assumption, the change; * tjS’i 

between elements are gradual (see the discussion m Sect. 5 of Ref. ^andtheerror 
introduced by a suitable stress averaging procedure should be small. In addition, St. 
Venant’s Principle is applicable to the model so that stress and strain states are only stn y 
valid in regions away from the blade ends and points where discontinuities in the blade 
geometry or structural properties occur. 
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The stress at such points may be approximated by either 1) referring the stress of 
each element to an ‘average reference frame’ and averaging these resulting stresses, or 2) 
merely averaging the stresses of each element. An average reference frame is here defined 
in terms of the anhedral and sweep that are averaged from the two adjacent elements. 
Presently, no averaging is done and the code simply chooses one of the elements and refers 
stresses in the local coordinates of that element 

The nodal deformations are determined from the modal amplitudes and any 
remaining steady state contributions from the blade rotation forces. Note that the actual 
displacement of each node also includes the rigid blade motion due to the cyclic pitch. 
However, this component of the motion does not contribute to the stresses since no 
deformations are associated with rigid body movement. Thus, from the original expression 
for the nodal displacements (in global coordinates), Equation 3.18: 

(q r ) =[X](T|) + {r 0c )e c + {q™) (3.40) 

where (r^) is the nodal displacement vector due to unit cyclic pitch and the residual: 


{q”*} = ([K*]- 1 - [XJtDrW) {F rot } (3.41) 

where, [D] = [ diag{ CO} 2 GM(i) ) ] is added to compensate for neglected modes in the 
analysis. This follows from the steady state deflection due to blade rotation obtained from 
the complete finite element analysis, 

{ q rot } = [K ir ]* 1 {F rot } (3.42) 

and the steady state deflection obtained by modal reconstruction. 


NMODE 


{q rot } = V {x(i)J — 

CO; 


1 


CO, GM(i) 


F R0T (i) 


i=l 

= [X] [diag{ * }] [X] T {F rot } (3.43) 

coj GM(i) 

As noted previously, this correction is necessary when the set of retained modes does not 
include an extension mode and thus is orthogonal to the rotational force vector. Without 
this adjusting term, the deflections due to the substantial rotational force would be 
unaccounted for thus resulting in substantial hub forces and axial stress errors. The actual 
deformations from Equation 3. 1 8 are 


(q) = M - {r c }0 c 

= [X]{Tl} + ( [Kjj ]' 1 - [XjtDrW ) (F 1-01 ) (3.44) 

These deformations are then referred to the local reference frames with Equations 3.38 and 
3.39 being used to obtain the blade stresses. 
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ProceduraUy, this is most easily accomplished in parallel with the construction of 
the array, SHAPES (nm,ir,k) which is the modal deflection, k, of mode, nm, at radial 
location, ir. Considering first the longitudinal stress: 


'XX 


E 1 


v. 


- [ TjcosP - £sin|3 ] { O 3 } 

- [ £cos|3 + T|sinJ3 ] {d > 3 } 

*P, X {*>2) 

i®2) 






([X] {r|} + {q 101 }) 


(3.45) 


The term in (•) is assumed to be referred to the local coordinate system of the element 
containing the point of interest. The full information necessary to determine c xx at any 
point on the cross- section is stored in the array SIGXX(nm,irJc), where. 

SIGXX(nm,ir,l) = - E {cos|3 {0 3 "} T : sinfl {<X> 3 } T } \ ■ [ 

[ x 8 Jnm 


SIGXX(nm,ir,2) = E { sinp {d> 3 ”) T : -cosp {0 3 "} T } \ : 

x 8 J 


nm 


x 9 


SIGXX(nm,ir,3) = E {0 2 ) T \ x io f SIGXX(nm,ir,4) = E {<X> 2 } 


T 


x ll 


' nm 


. t r x 1 2 1 

SIGXX(nm,ir,5) = E {<D 2 } M x 13 ( 

l x 14 Jnm 


x 9 ] 
x 10 \ 
x 1 1 J nm 

(3.46) 


where (x^-.x^} 1 are the local coordinate deflections due to mode, nm. In particular, if 

the global deflections corresponding to mode, nm, for element, ie, are, {Xj.. .X14 } (these 
are simply those D.O.F. of the modal vector, (X), that correspond to element, ie) then, 
{xi ...x 14 } t = [T rot ] T {X 1 ...X 14 ) T . Then at any cross-section point, (n,C), at radial 
location, ir, the longitudinal stress generated by mode, nm, is 

o M (nm) = -q SIGXX(nm,ir,l) + C SIGXX(nm,ir,2) 

+ 'P.x SIGXX(nm,ir,3) + 'P SIGXX(nm,ir,4) 

+ SIGXX(nm,ir,5) (3-47) 

and the total stress, 

NMODE 

c xx = (^(0) + X °xx( nm ) (3 48) 

nm=l 

where the constant stress: 

a xx (0) = T| SIGXX(0,ir,l) + C SIGXX(0,ir,2) 

+ 'P.x SIGXX(0,ir,3) + *P SIGXX(0,ir,4) 

+ SIGXX(0,ir,5) ( 3 - 49 ) 

and, 
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' -rot ' 


SIGXX(0,ir,l) = - E {cosp {0 3 "} T : sinp {0 3 "} T } \ 




SIGXX(0,ir,2) = E { sinp {<D 3 "} T : -cosP {<D 3 "} T } 


r q7° 


SIGXX(0,ir,3) = E {<X> 2 ’} T S 


SIGXX(0,ir,4) = E {<D 2 "}M 
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q 8 


(3.50) 


A similar procedure may be followed for the shear stresses, a XT ^ and o x £ . However, 
shear stresses are not currently implemented in RotorCRAFT. ^ 

In each case, the evaluation of the stress requires the specification of a particular 
point on the cross-section of the rotor blade. Typical strain gauge installations involve 
differencing the measured strain from different locations on die blade. Similarly, blade 
stress is most directly computed by evaluating the local stress at different locations on the 
blade and differencing the resulting predictions. Examples of such calculations are 
discussed in Section 5. 
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4. OVERVIEW OF CODE OPERATION 

This section provides a brief summary of the trim procedures embedded in 
RotorCRAFT (Mod 1.0), with particular reference to the new provisions for arbitrary 
neriodic Ditch input Also, the discussion below provides an outline of the preliminary 
interface ^^eloped here for interaction of RotorCRAFT with NASA's WOPWOP noise 
prediction code. It is emphasized at the outset that the present interface produces results 
that are considered suitable only for demonstration calculations and preliminary software 
development, not for correlation of acoustic data. The development of a more refined and 
complete analysis of rotor noise is the subject of ongoing research, the preliminary results 
of which are reported in References 25 and 26. 

4.1 Rotor Trim Procedures in the Presence of Higher Harmonic Pitch 

As discussed above and in Reference 2, the Mod 0.0 variant of the RotorCRAFT 
analysis was designed to solve for the aerodynamic loads on isolated helicopter rotors in 
steady forward flight. The rotor hub is assumed to be fixed in space and the rotor shaft is 
assumed to be oriented at a fixed angle with respect to the free stream. The solution 
method is directed at obtaining periodic solutions for the wake geometry, the blade motion 
the aerodynamics loads, and now, in the Mod 1.0 variant, for the blade stresses and hub 
loading. This objective was judged to be appropriate for this particular effort, given the 
focus on steady-state forward flight. The "transient" solution achieved during the iteration 
process does not represent a time-accurate calculation of rotor loading; only the converged 
result can be considered physically valid and consistent. The approach used to obtain this 
converged result is now described. 

The outermost loop in the iteration process is the computation of the evolution of 
the wake geometry over one rotor revolution with blade motion and aerodynamic loading 
fixed The simulation is ordinarily run for a discrete number of rotor revolutions, and 
between each revolution the cyclic pitch inputs and the collective are adjusted to meet the 
conditions specified by the user for thrust and first harmonic flapping response (ordinarily, 
zero first harmonic flapping is desired, though many experimental cases in the literature 
have been run with nonzero flapping). This trim procedure represents the second loop 
within the overall iteration. Finally, the innermost loop is the calculation of the blade 
motion that is consistent with the current values of the pitch control settings; this calculation 
is carried out using the dynamic model described in Section 3 above and is performed 
assuming the wake-induced downwash is fixed at the most recent estimate obtained from 
the outer loop involving the wake evolution. 

The principal change in this procedure in the Mod 1.0 variant is that the pitch 
command may now be of the form 

0(y) = 0 O + 0 lc cos y + 0i s sin y + ...+ 0 n cCos n y + 0nsSin ny + ... (4.1) 

The root pitch may thus be adjusted to any desired set of harmonic inputs. In addition, the 
user may also select an arbitrary time history of root pitch inputs with the appropriate 
parameter selection, with the constraint that it be periodic with the rotor rotation. If the user 
wishes, the trim procedure can be switched off so that the cyclic pitch time history is fixed 

according to the specified parameters 0no 6ns- If not, the trim procedure will still leave the 

higher harmonic content of the cyclic pitch input unchanged, (0 n c and e ns for n > 1), but 

will modify 0i c and 0i s until the user requested first harmonic rigid flapping values are 

obtained. 
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The first step in starting the overall solution procedure is to determine an initial 
estimate for the blade motion and the cyclic controls using the vortex lattice method and an 
estimate of wake-induced velocity based on a simple, prescribed downwash distribution. 
This calculation includes the higher harmonic pitch input, but can be bypassed in favor of 
an arbitrary periodic input, as noted above. Using the blade motion and bound circulation 
estimates so generated, the first loop of the analysis proceeds, beginning with a kinematic 
wake whose geometry is determined by the free stream and the downwash at the rotor disk. 
The wake then is allowed to evolve for one rotor revolution, during which the time history 
of the wake-induced flow field at each evaluation point on the rotor blade is computed and 
stored. At the end of this revolution, this time history of wake-induced flow is passed to a 
trim and dynamics routine that updates the control settings and the blade motion solution to 
be consistent with the free wake flow field. 

If automatic trim has been selected, the program uses the updated velocity field and 
checks whether the first harmonic components of rigid flapping are within a specified 
tolerance of the desired level. If not, the cyclic pitch is adjusted and the blade motion 
calculation repeated until they are. Once trim has been achieved, the thrust is checked to 
see if it lies within a specified tolerance of the desired level. If not, the collective is 
adjusted and the steps just above are repeated until this condition is met. The iteration 
history for thrust coefficient is tracked and previous results are used to accelerate the 
convergence process. This process has proved to be sufficiently robust to handle a wide 
variety of rotor configurations, as will be evident from the sample cases presented in both 
Reference 2 and in Section 5 below. A flow chart of the sequence of events in the trim 
cycle is given in Figure 4-1; more information in the input requirements can be found in 
Reference 15, and a detailed discussion of the overall trim calculation can be found in 
Reference 2. 

Once the blade motion solution has converged and the final aerodynamic loading 
distribution has been determined, the information needed for the computation of rotor blade 
hub loads and internal stresses is available. The computation of the hub loads proceeds as 
described in the previous section and produces time histories of the forces and moments 
exerted on the rotor hub, with the quantities expressed in the global (shaft) reference frame. 
The computation of the internal blade stress requires the specification of the radial and 
vertical positions of evaluation points relative to the blade's elastic axis (Ref. 15). The 
calculations described in Section 5 will include representative results of both stress and hub 
load computations. 

4.2 WOPWOP Coupling 

The development of comprehensive rotor noise prediction capabilities is a topic of 
ongoing research, and the development of efficient flow field prediction methods to support 
this work was the topic of recent NASA-sponsored research (Ref. 26). As a first step in 
adapting the flow field prediction code described in Reference 26 to the problem of noise 
prediction, a preliminary scheme has been developed for generating the unsteady blade 
surface pressure information required to carry out noise computations. 

The primary objective of the present limited effort was to produce output from 
RotorCRAFT compatible with the input requirements of an independent acoustic program, 
WOPWOP, developed at NASA Langley (Ref. 27). WOPWOP was selected because of its 
robustness in analyzing general helicopter rotor configurations. The fundamentals of its 
acoustic analysis are based on Farassat's subsonic formulation of the Fowcs Williams- 
Hawkings equation derived for the noise generated by surfaces in arbitrary motion. 
Assuming that the surface loading distribution is known beforehand, Farassat's scheme 
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provides an effective prediction of noise generated by surfaces in subsonic flow. The 
scheme implemented in WOPWOP was specially tailored for the helicopter rotor acoustic 
problem, and both the blade motions and surface pressures are required inputs to the 
program. The surface pressures, in particular, must be accurately specified to a high 
degree of resolution to obtain good noise predictions. 

Since the Mod 0.0 version of RotorCRAFT focused largely on the prediction of 
aerodynamic loading that contributes to vibratory airloads, a substantial increase in 
refinement is required to determine a detailed distribution of chondwise airloads on the rotor 
blades. Substantial research work to date has been dedicated to develop practical models 
for die accurate solution of high-frequency surface airloads distribution for acoustics 
applications. Since many issues pertaining to the creation of a generally applicable model 
for unsteady loading remain unresolved, it was judged appropriate for current purposes to 
adopt a relatively simple model based on the selection of an aerodynamic transfer function 
to predict the lift response to the unsteady upwash field experienced by the rotor blades. 
The approach taken here is very similar to that described in Reference 26. 

A routine designed to provide the necessary interface has been implemented in 
RotorCRAFT (Mod 1.0) as a post-processor to the rotor- wake calculation and is invoked 
only after both the blade dynamics and the flow field have been converged (Fig. 4-2). The 
converged unsteady upwash predicted on each of rotor blades (composed of the time- 
varying free stream, the companion blades and the unsteady wake) is then treated as an 
arbitrary gust. To find the lift response, each blade is segmented in the spanwise direction 
and each segment is treated as a two dimensional flat plate airfoil. The Kussner function is 
used with the Duhamel superposition integral to evaluate the loading response of the airfoil 
to this simplified gust problem. The basic transfer function approach assumes a two 
dimensional linear wake extending far downstream of the airfoil. However, to approximate 
the effect of the skewed helical wake of a rotor blade, which folds back on itself rather than 
extending to infinity, adjustments were made to the analysis being used that effectively 
truncated the wake at a finite distance downstream of the blade. Both the amplitude of the 
pressure response and the nondimensional time were corrected for compressibility in the 
spirit of Prandtl-Glauert. This approach clearly requires many approximations, but it does 
provide a consistent estimate of the unsteady surface pressure distribution on the rotor 
blades suitable for direct input to WOPWOP. 

This method for determining the surface airloads distribution on the rotor blade is 
sufficient for demonstrating the feasibility of the transfer function approach to helicopter 
rotor aeroacoustic analysis. As it was discussed in Reference 26, calculations thus 
performed can provide some physical insight into the qualitative sensitivity of the predicted 
loading noise to selective input parameters and so help define some issues of interest for 
follow-on efforts. However, given the many approximations inherent in the present 
approach, it is not considered appropriate for the correlation of the acoustic predictions 
quantitatively with experimental data. It is recognized that substantial advances in the 
evaluation of the distributed surface loadings of the rotor blades have to be achieved before 
such correlation can be made with hope of success. Among the important issues to be 
addressed in improving this mode is the fact that the current analysis is based on an inviscid 
calculation and viscous shear stresses experienced on the blade surface have been 
neglected in the airloads evaluation. In addition, effects of three-dimensionality are in 
general of considerable importance in computing surface airloads, especially at the blade 
tip, and thus a transfer function with explicit representation of finite aspect ratio effects 
must be developed. Finally, efficient prediction of wake-induced inflow with high 
temporal resolution is also required. Such issues define some of the major challenges to be 
faced in converting the present analysis to a comprehensive analysis of rotor acoustics. 
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Figure 4-2. Implementation of the preliminary interface between 
RotorCRAFT and WOPWOP for the analysis of rotor noise 
predictions. 
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5. DATA CORRELATION STUDIES 

Several data correlation studies using the baseline RotorCRAFT code were carried out 
as part of the development of the Mod 0.0 version of the code (Ref. 2). These calculations 
focused almost exclusively on steady and unsteady aerodynamic loading, though limited 
correlation of wind tunnel data on rotor flapping at low speed was also conducted. For the 
present effort, blade stress and hub load computations were of primary interest 

5.1 Blade Stress Calculations 


. . _P” e convenient and widely used source of blade stress data is the data base from the 
H-34 flight tests described in Reference 28. The correlation studies undertaken here will 
focus on examining flight conditions from that test. Undertaking blade stress calculations 
first requires die selection of mass and stiffness properties for the blade. Reference 28 did 
not contain this information, but it was assumed that the documentation of the H-34 blades 
“ d mnnel test in Reference 29 could be used. The blade structural calculations 

m RotorCRAFT require that this information be input in the form of radial distributions of the 
following quantities: Young's modulus and shear modulus; cross-sectional area; offset of the 
cross-sectional area centroid from the elastic axis; second moment of area about elastic axis - 
i>t. Venant s torsion constant for the section; angle between the local elastic axis and global 
span wise axis; torsion/bending warping coefficients; mass per unit length; sectional center of 
mass offsets from the elastic axis; moments of inertia per unit length; and the cross product of 
inertia at each section. Detailed discussion of the input requirements for the RotorCRAFT 
structural properties file are given in Reference 15. 


d , th , e ™P ortant considerations in the analysis of blade stress using the 

otorCRAFT analysis (or any comparable finite element analysis) is appropriate selection of 
the number and sizing of the elements used to discretize the blade. As noted above the 
structural model consists of beam elements with fourteen degrees of freedom. Numerical 
experimentation with the H-34 configuration indicated that ten elements were sufficient to 
resolve the stnictural deformation for the purposes of blade stress calculations. More refined 
input files with as many as nineteen elements were prepared and used as input, but the results 
were found to be negligibly different from the ten-element cases discussed below. 

The stress calculation procedures described in Section 3 above require the selection of 
a point of evaluation for the bending stress, described by a radial position and a vertical 
displacement from the blade’s elastic axis. In the beam-rod approximation used in the 
present structunti model, this corresponds to finding the stress in particular fibers at specified 
distances away from the elastic axis. The net stress experienced at a particular radial station 
can be found by differencing the stress computed at measurement points placed above and 
below the elastic axis; this differencing serves to subtract out the large mean stress 
component due to rotational effects, leaving an unsteady component attributable to 
aerodynamic loading and structural deflection. 


The primary case of interest here is the advance ratio 0.29 case from Reference 28. 
Figure 5-1 shows a comparison of the total unsteady bending stress computed at r/R=0 65 
with the measured value. There is some phase error in the first two quadrants, but the 
agreement is on the whole quite close. Of more general interest, however, is the vibratory 
blade stress, i e. the components at 3P and higher for this four-bladed rotor. Recent work by 
Bousman (Ref. 30) involved general characterization of blade stress measurements across a 
wide range of rotor configurations, including the H-34 test examined here. The results 
quota! there indicated that the flapwise vibratory stress response of a variety of blades 
including the H-34 was dominated by a 3P component. Figure 5-2 shows this behavior for 
several different rotors, while Figure 5-3 features a surface plot of such response, drawn 
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Figure 5-1. Total unsteady (harmonics 1-10) blade bending stress for the H-34 at 
advance ratio 0.29 . 
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Figure 5-2. General trends in bending stress data for four rotors in high speed 
flight (from Ref. 28). 
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Figure 5-3. Surface plot of harmonics 3-5 of measured blade bending stress for the 
H-34 at advance ratio 0.29 (from Ref. 26). 
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Figure 5-4. Vibratory (harmonics 3-10) blade bending stress for the H-34 at 
advance ratio 0.29 : r/R = 0.65 . 







from the data in Reference 28. Figures 5-4 to 5-6 show the predicted levels of stress for 
harmonic components from 3P to 10P. All the results show good qualitanve c^elaaons 
with the levels of response being particularly well predicted at the stetions r/R^575 and 
0.65. The results farther inboard tend to underpredict the level of stress response, 
particularly in the first quadrant. 

5.2 Hub Load Calculations 

As noted in Section 1, the experimental results reported in Reference 1 indicate 
considerable promise in the application of higher harmonic pitch control for nm^ r^uchom 
Several of the test cases described in Reference 1 were examined using RotorCRAFT (Mod 
1.0) to assess the correlation between the predicted S 

rotor was tested in the Transonic Dynamics Tunnel (TDT) at ^SM*nriey _and ^vnll be 
referred to hereafter as the TDT rotor. The rotor radius was 4.58 ft. and its chord was 
constant over the span at a value of 0.35 ft. The bl^es were untwisted, and ^ 
sectional properties were supplied by NASA personnel and suitably arranged for inputmto 
RotorCRAFT (Mod 1.0). The rotor was tested in heavy gas (freon) so that Mach-scale tests 
could be conducted at reduced tip speeds. The tip speed used in the tests m Reference 1 was 
31 1 fps, corresponding to a rotor rotation frequency of 68 rad/sec. The tests to be examined 
in detail here are four flight conditions from those documented m Reference 1, these 
particular cases were chosen since they ’bracket’ the advance ratio and shaft angle range 
considered. The four flight conditions are: advance ratio 0.166 at shaft angle of attack 
+6.0°, (designated 'Case A', corresponding to Figure 17a of Reference 1); advance ratio 
0.30 at shaft angle +2.0° (Case C); advance ratio 0.166 at shaft angle 0 (Case G); and 
advance ratio 0.3 at shaft angle -4.0° (Case I). In each case the rotor thrust coefficient was 

0.0050. 

The first step in this correlation study was to compute the modal properties of the 
blade using the finite element structural analysis in the RotorCRAFT code, as described in 
Reference 2. The predicted modal frequencies for the first three out-of-plane tendrng modes 
were 1.05P, 2.62P, and 4.90P, while in-plane frequencies were computed to be 03 IP (rigid 
lag) and 4.12P (first elastic lag). The corresponding frequencies predicted by CAMKAD 
were 1.04P, 2.61P, 4.82P, 0.29P and 4.1 8P, respectively (Ref. 31). The modjd frequency 
calculations also predicted that the first elastic torsion frequency is quite high (8.60P). Since 
the primary experimental data of interest is 4P hub loading, it appears ^likely that torsional 
deflection participates to any large degree in this response. RotoiCRAFT (Mod 1.0) permits 
the user to select which elastic deflection modes are to be included in the analysis. It was 
determined that an appropriate numerical model featured three out-of-plane bending modes 
(rigid flapping and the first two elastic bending modes) as well as two m-plane modes (rigid 

lag and first elastic lag). 

In addition to the specification of the blade’s geometry and structural properties, the 
RotorCRAFT code requires the user to select a vortex lattice configuration for the blade as 
well as a wake model suitable to the flight condition. For all the calculations described 
below the blade model uses a vortex lattice with 30 vortex quadrilaterals spanwise and one 
chordwise quad. For Case C and Case I, the calculations use a maximum of eighteen wake 
filaments trailing from Zone 1 and ten from Zone 2 (see Ref. 15 for definitions of these 
Zones), with one full turn of full-span CVC wake and one turn of free wake extensions m 
each case. For the relatively low speed calculations of Case A and Case G, a longer free 
wake was required, though the same number of filaments were trailed from the span. In each 
case, default selections for the vortex cores were made. It is possible that changes of one or 
more of these parameters may affect the results obtained below, though these selections were 
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Figure 5-5. Vibratory (harmonics 3-10) blade bending stress for the H-34 at 
advance ratio 0.29 : r/R = 0.57 . 



Figure 5-6. Vibratory (harmonics 3-10) blade bending stress for the H-34 at 
advance ratio 0.29 : r/R = 0.45 . 




made on the basis of substantial numerical experimentation with other similar rotor 
configurations. 


Figure 5-7 shows a snapshot of wake geometry in Case I both before and after the 
application of 4P higher harmonic pitch. The change in wake structure is quite marked; the 
presence of closed circles in the wake indicate the appearance and disappearance of loading 
peaks due to the time-varying pitch. Since the flow field produced by the is 

correspondingly altered, this illustrates one of the advantages of using the refined CVC wake 
model for this type of calculation. 


Reference 1 specified the higher harmonic pitch input used in each case. As noted 
above, RotorCRAFT (Mod 1.0) can take as input arbitrary periodic pitch time histones 
including higher harmonic forms such as : 


0(y) = e 0 + OicCos y + OisSin y + ...+ 0 nc cos ny + e ns sin ny + 


(5.1) 


For these correlation calculations, the steady and first harmonic pitch inputs were chosen to 
trim the rotor to zero first harmonic flapping at the desired thrust coefficient. All other 
harmonic components were set to zero, except the 4P terms. Reference 1 specified these 

inputs in terms of a magnitude 0 C and a phase angle y c . These were transformed into the 
required components of 0 using 


0 4c = 9 C cos 4y c and 0 4s = 0 c sin4y c 


These inputs replicate the 4P collective pitch inputs used in the test, which featured 
amplitudes of 0 C = -0.6° and 0 C = -1.2°. The value of the phase angle y c ranged from 0° to 
90° for each of the flight conditions examined. 


Before calculating the vibratory load response to higher harmonic control inputs 
(HHC) the prediction of levels of vertical hub force with no HHC was first considered. The 
data presented in Figure 5-8 were drawn from Figure 17 of Reference 1 -Except for Case C, 
the load magnitudes predicted by RotorCRAFT (Mod 1.0) without HHC agree with the 
measured data to within 10%. Therefore, any discrepancies greater than this that are 
encountered in the comparisons that follow are assumed to be associated with the effect of 
modeling higher harmonic control inputs. 


The tests described in Reference 1 were performed on a non-scaled rotor hub that 
generated substantial inertial vibratory load when subjected to the 4P higher harmonic pitch 
input. The quantities measured in the test were the forces and moments transmitted to the 
stationary reference frame. For comparison with RotorCRAFT (Mod 1.0) predictions, the 
inertial forces associated with the swashplate motion were estimated and subtracted from the 
total forces obtained during the test. The aerodynamic and aerodasac force contnbutions 
remaining are equivalent to those calculated by the RotorCRAFT (Mod 1.0) code. The hub 
force data presented in Figures 5-9 and 5-10 were obtained from the following formula 


(F 4 p)aero = (E 4 p)total ’ F^pIhHC 

where Fi is the magnitude of the inertial force associated with the swashplate motion. The 
vector sense of the forces here represents the sin 4y and cos 4y components; the inertial 
force is in phase with the pitch control acceleration. A value of Fi = 44.74 lbs was used for 
the current analysis (Ref. 32). 
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a) With no HHC input 



b) With HHC input 


Figure 5-7. Wake geometry for ji = 0.3, a s = -4.0° with no HHC pitch input (top) 

and with HHC pitch input (0 C = -1.2% V c = OP) (bottom). (Vertical 
scale expanded by a factor of 5.0 .) 
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